#First we load our data into our environment

library(tidytuesdayR) 

tuesdata <- tidytuesdayR::tt_load('2025-06-24') 

cases_month <- tuesdata$cases_month #Loading datasets into environment
cases_year <- tuesdata$cases_year   #Loading datasets into environment
#We start by loading our packages

library(rayshader) 
library(rgl)
library(av)
library(gifski)
library(magick) #The first five are just for the 3d-model gif. 
library(ggplot2)
library(dplyr)
library(magrittr) #should be included in dplyr, but didnt work?
library(plotly) #For interactive ggplot
#Looking at variable names
names(cases_month)
##  [1] "region"                "country"               "iso3"                 
##  [4] "year"                  "month"                 "measles_suspect"      
##  [7] "measles_clinical"      "measles_epi_linked"    "measles_lab_confirmed"
## [10] "measles_total"         "rubella_clinical"      "rubella_epi_linked"   
## [13] "rubella_lab_confirmed" "rubella_total"         "discarded"
names(cases_year)
##  [1] "region"                                                         
##  [2] "country"                                                        
##  [3] "iso3"                                                           
##  [4] "year"                                                           
##  [5] "total_population"                                               
##  [6] "annualized_population_most_recent_year_only"                    
##  [7] "total_suspected_measles_rubella_cases"                          
##  [8] "measles_total"                                                  
##  [9] "measles_lab_confirmed"                                          
## [10] "measles_epi_linked"                                             
## [11] "measles_clinical"                                               
## [12] "measles_incidence_rate_per_1000000_total_population"            
## [13] "rubella_total"                                                  
## [14] "rubella_lab_confirmed"                                          
## [15] "rubella_epi_linked"                                             
## [16] "rubella_clinical"                                               
## [17] "rubella_incidence_rate_per_1000000_total_population"            
## [18] "discarded_cases"                                                
## [19] "discarded_non_measles_rubella_cases_per_100000_total_population"
#Recoding regions for OUR understanding
cases_year <- cases_year %>% 
  mutate(region = recode(region, `AFRO` = "African Region", `AMRO` = "Regon of Americas", `EMRO` = "Eastern Mediterraenean Region", `EURO` = "European Region", `SEARO` = "Sout-East Asian Region", `WPRO` = "Western Pacific Region"))
cases_month <- cases_month %>% 
  mutate(region = recode(region, `AFRO` = "African Region", `AMRO` = "Regon of Americas", `EMRO` = "Eastern Mediterraenean Region", `EURO` = "European Region", `SEARO` = "Sout-East Asian Region", `WPRO` = "Western Pacific Region"))
#First we want to look at the overall average incidence rate pr. year worlwide
cases_year %>%
  group_by(year) %>% #Grouping by year
  mutate(avg_reg = mean(measles_incidence_rate_per_1000000_total_population)) %>% #Creating new variable with mutate
ggplot(aes(x = year, y = avg_reg)) + #inserting mapping-argument in overall ggplot for less typing
  geom_point() +
  theme_bw() +
  labs(title = "Avg. measles incidence rate worldwide \n(pr. 1,000,000)", x = "Year", y = "Avg. Incidence Rate") + #Notice the \n for changing                                                                                                                      #line
  geom_smooth()

#Now we want to facet by region. We supress warnings and messages in the output in the first row. 
plotfacet <- cases_year %>%
  group_by(year, region) %>%    #We group by year AND region so we get an annual avg. pr. region
  mutate(avg_reg = mean(measles_incidence_rate_per_1000000_total_population)) %>% 
ggplot() + 
  geom_point(aes(x = year, y = avg_reg, colour = region, shape = region)) +
  theme_bw() +
  labs(title = "Avg. measles incidence rate pr. region pr. year(pr. 1,000,000)", x = "Year", y = "Avg. Incidence Rate", color = "Region", shape = "Region") +
  geom_smooth(aes(x = year, y = avg_reg), color = "grey", se=FALSE, show.legend = FALSE) + 
  facet_wrap(~ region, scales = "free") +
  geom_line(aes(x = year, y = avg_reg, colour = region, shape = region))

ggplotly(plotfacet)
cases_year1 <- cases_year %>% 
  group_by(year, region) %>% 
  mutate(avg_reg = mean(measles_incidence_rate_per_1000000_total_population)) 

hexgg <-  ggplot(cases_year1, aes(x =  year, y = avg_reg)) +
    stat_bin_hex(aes(fill = after_stat(density), colour = after_stat(density)), 
                 bins = 10,
                 linewidth = 1) +
    scale_fill_viridis_c(option = "B") +
    scale_color_viridis_c(option = "B", guide = "none") +
    labs(x = "Year", y = "Avg. incidence rate pr. region", fill = "",
         colour = "") +
    theme_minimal()
hexgg

cases_month %>%
  mutate(region = as.factor(region)) %>% 
  ggplot(aes(x = measles_suspect, y = measles_lab_confirmed, color = region)) +
  geom_point() + 
  facet_wrap(~region) + 
  theme_bw() + 
  labs(x = "Suspected cases of measles", y = "Laboratory confirmed cases of measles", color = "Region")